Highly accurate estimates of the high-spin/low-spin energy difference ΔE_HL^el in the high-spin complexes [Fe(NCH)₆]²⁺ and [Co(NCH)₆]²⁺ have been obtained from the results of CCSD(T) calculations extrapolated to the complete basis set limit. These estimates are shown to be strongly influenced by scalar relativistic effects. They have been used to assess the performances of the CASPT2 method and of 30 density functionals of the GGA, meta-GGA, global hybrid, RSH and double-hybrid types. For the CASPT2 method, the results of the assessment support the proposal [Kepenekian, M.; Robert, V.; Le Guennic, B. J. Chem. Phys.2009, 131, 114702] that the ionization potential–electron affinity (IPEA) shift defining the zeroth-order Hamiltonian be raised from its standard value of 0.25 au to 0.50–0.70 au for the determination of ΔE_HL^el in Fe(II) complexes with a [FeN₆] core. At the DFT level, some of the assessed functionals proved to perform within chemical accuracy (±350 cm^-1) for the spin-state energetics of [Fe(NCH)₆]²⁺, others for that of [Co(NCH)₆]²⁺, but none of them simultaneously for both complexes. As demonstrated through a reparametrization of the CAM-PBE0 range-separated hybrid, which led to a functional that performs within chemical accuracy for the spin-state energetics of both complexes, performing density functionals of broad applicability may be devised by including in their training sets highly accurate data like those reported here for [Fe(NCH)₆]²⁺ and [Co(NCH)₆]²⁺.